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A WWW interface for the simulation of spectral energy distributions of optically 
^ ■ thin dust configurations with an embedded radiative source is presented. The den- 

sity distribution, radiative source, and dust parameters can be selected either from 
an internal database or defined by the user. This tool is optimized for studying cir- 
cumstellar debris disks where large grains (a gra i n 3> l^m) are expected to determine 
the far-infrared through millimeter dust reemission spectral energy distribution. The 
tool is available at http://aida28.mpia-hd. mpg.de/~swolf/dds. 

Key words: Radiative transfer, dust, absorption, scattering, debris disk - 
astrophysics 



Preprint submitted to Elsevier Preprint 



2 February 2008 



PROGRAM SUMMARY 



Title of program: Debris Disk Radiative Transfer Simulator (DDS) 
Catalogue number: (supplied by Elsevier) 
Program obtainable from: 

CPC Program Library, Queen's University of Belfast, N. Ireland (see application 
form in this issue) 

Licensing provisions: none 

Computers: 

PC with Intel(Fi) XEON(TM) 2.80 GHz processor 

Operating systems under which the program has been tested: 
SUSE Linux 9.1 

Programming language used: 

Fortran 90 (for the main program; furthermore Perl, CGI and HTML) 

Memory required to execute with typical data: 10 s words 
No. of bits in a word: 8 

No. of lines in distributed program, including test data, etc.: 
9264 

Keywords: 

Nature of physical problem 

Simulation of scattered light and thermal reemission in arbitrary optically dust 
distributions with spherical, homogeneous grains where the dust parameters (op- 
tical properties, sublimation temperature, grain size) and SED of the illuminat- 
ing/heating radiative source can be arbitrarily defined (example application: Wolf 
& Hillenbrand 2003). The program is optimized for studying circumstellar debris 
disks where large grains (i.e., with large size parameters) are expected to determine 
the far-infrared through millimeter dust reemission spectral energy distribution. 

Method of solution 

Calculation of the dust temperature distribution and dust reemission and scattering 
spectrum in the optically thin limit. 

Restrictions on the complexity of the problem 

1) The approach to calculate dust temperatures and dust reemission spectra is only 
valid in the optically thin regime. The validity of this constraint is verified for each 
model during the runtime of the code. 

2) The relative abundances of different grains can be arbitrarily chosen, but must 
be constant outside the dust sublimation region., i.e., the shape of the (arbitrary) 
radial dust density distribution outside the dust sublimation region is the same for 
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all grain sizes and chemistries. 

3) The size of upload files (such as the dust density distribution, optical data of the 

dust grains, stellar spectral energy distribution, etc.) is limited (see http:/ /aida28. mpia-hd.mpg.de/~swolf/dds/ 
for current file size limits). However, the resulting limitation to the complexity of 
possible model definitions is marginal only. 

Typical running time 

3 sec - 30 min (depending on the complexity of the model) 
Unusual features of the program 

The program as provided through the CPC Program Library is equipped with an 
HTML user interface. It is installed and available at 
http : / / aida28 . mpia-hd . mpg . de/ ~swolf /dds 

References 

[1] Wolf, S., Hillenbrand, L.A., 2003, Astroph. Journal, 596, 603 
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LONG WRITE-UP 



1 Introduction 

The Debris Disk Radiative Transfer Simulation Tool (DDS) was developed 
with the aim to provide a flexible tool for the simulation of spectral en- 
ergy distributions (SEDs) of circumstellar debris disks. Debris disks are solar 
system-sized dust disks with micron-sized grains produced as by-products of 
collisions between asteroid-like bodies being left over from the planet forma- 
tion process. In the case of our solar system, the debris of Jupiter-family 
short-period comets and colliding asteroids represents the dominant sources 
of zodiacal dust located between Mars and Jupiter. A second belt of dust is 
located beyond the orbit of Neptune (see, e.g., Dermott et al. 1992, Liou, Der- 
mott, & Xu 1995). Besides the solar system, optical to mid-infrared images of 
(3 Pic (see, e.g., Kalas & Jewitt 1995; Weinberger, Becklin, & Zuckerman 2003) 
and submillimeter images of Vega, Fomalhaut, and e Eri (Holland et al. 1998; 
Greaves et al. 1998) have revealed spatially resolved debris disks which were 
first inferred from observations of infrared flux excesses above photospheric 
values with IRAS. The mass of small grains in debris disks and therefore 
the thermal dust reemission from these disks is expected to be much smaller 
than in the case of TTauri disks. For this reason, only a very limited sam- 
ple of observations exists so far. However, because of the high sensitivity of 
the mid-infrared detectors aboard the Spitzer Space Telescope followed by the 
Stratospheric Observatory for Infrared Astronomy (SOFIA), a substantial in- 
crease in the total number and in the specific information about debris disks 
is expected (c.f. Meyer 2002). 

The range of applications of the DDS, however, extends far beyond debris 
disks (or similar, e.g. disk-like structures) allowing to derive the reemission 
and scattered light SED for any three-dimensional optically thin dust distri- 
bution: In the optically thin limit only the radial distance r between each 
individual grain of a dust distribution determines this grain's contribution 
to the net SED. Thus, any arbitrarily shaped configuration described, e.g. 
by a three-dimensional denisty distribution p(x, y, z) can be reduced to an 
equivalent one-dimensional distribution p(r) in respect of its resulting ther- 
mal reemission and scattered light SED (r is the radial distance from the 
central star / heating source). The DDS is therefore best characterized as a 
tool for the simulation of scattered light and thermal reemission in arbitrary 
optically thin dust distributions with spherical, homogeneous grains where the 
dust parameters (optical properties, sublimation temperature, grain size) and 
SED of the illuminating/heating radiative source can be arbitrarily defined. 
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2 The method 



The solution of the radiative transfer problem implemented in the DDS is 
based on three assumptions: 

(1) The dust configuration is optically thin along any line of sight for both 
the stellar as well as for the reemission radiation from the dust, 

(2) The dust grains are compact spheres with a homogeneous chemical struc- 
ture, and 

(3) The dust grains are large enough to be in thermal equilibrium with the 
ambient radiation field. 

In the optically thin limit, each dust grain is heated by direct stellar radiation 
only Thus, the dust grain temperature is a function of the optical parame- 
ters of the grains, the incident stellar radiation, and the distance r from the 
star. In this case the radiative transfer equation has a simple solution which 
allows one to derive the distance from the star at which the dust has a certain 
temperature. Let 

L A , S = 4nR 2 s nB x (T s ) (1) 



be the monochromatic luminosity of the star (radius R s , effective temperature 
T s ) at wavelength A and 

L t = L ^ bs ^ and ( 2 ) 

= Ana'Qf s nB x (T g ) (3) 



be the absorbed and reemitted luminosity of a dust grain with radius a and 
resulting temperature T g at the (unknown) radial distance r from the star. 
Using the constraint of energy conservation 

oo oo 

/ L^dX = J L^dX (4) 





one derives the distance of the grain from the star as 



r(T.) = £. 



SS°Qf(a)B x (T a )dX 



Q*»(a)B x (T g )d\' 



(5) 



If the dust sublimation temperature is known, Eq. 5 also allows one to estimate 
the sublimation radius r su b for each dust component in the shell (characterized 
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by the grain radius and chemical composition). The flux of the light scattered 
by a single dust grain amounts to 



where A is the dust grain's albedo. The net spectral energy distribution results 
from a simple summation of the reemitted and scattered light contributions 
from all grains. 

The interaction of the stellar radiation field with the dust grains - charac- 
terized by the efficiency factors Qf> s and (and thus the albedo A) - is 
described by Mie scattering theory in the DDS. The Mie scattering function is 
calculated using the numerical solution for the estimation of the Mie scatter- 
ing coefficients published by Wolf & Voshchinnikov (2004; see also Voshchin- 
nikov 2004), which achieves accurate results both in the small as well as in 
the - arbitrarily - large size parameter regime. 




(6) 
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3 Numerical Implementation 



® WWW interface 

> General model description 

> Possibility to upload 

- SED of radiative source 

- Spatial dust density distribution 

- Dust data 




© OTC 

(Radiative transfer tool for 
optically thin configurations) 



©WWW output 

> Resulting SED 

> Simulation by-products 
(Dust scattering/absorption/etc. efficiency factors, 
Dust albedo, Radial temperature distribution, etc.) 

> Simulation protocol 

® WWW interface 

> Possibility to modify selected model parameters 



Fig. 1. DDS data flowchart. 

In this section the general scheme and main features of the DDS are described. 
For this reason a simplified chart symbolizing the data flow is given in Fig. 1. 
The single units, symbolized by letters A - F in Fig. 1, solve following tasks: 

A: 

WWW interface: Model setup (see Fig. 2). 

Each debris disk model to be investigated with the DDS (or dust configuration 
model in general) is defined by 

(1) The SED of the radiative source which is given 

(a) analytically as a blackbody radiator characterized by its temperature 
and total luminosity, 

(b) by an internally defined SED (e.g. the solar SED), or 

(c) by the user; the format of the upload-files containing the SED is given 
in Sect. A.l. 



- © Internal Database 

- SEDs for radiative source 

- Dust data 
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Image not available 



Fig. 2. DDS: Input mask. Figure available in the complete article - see: 
http: / / aida28.mpia-hd.mpg.de/~swolf/dds / doc / swolf_dds.pdf 

(2) The inner and outer radius of the disk. If an uploaded density distribution 
(see 3. below) is considered, the inner and outer radius defined in that 
file are used. If a fixed inner radius is chosen which turns out to be within 
the sublimation region of a particular dust species (defined by its grain 
size and chemical composition), the sublimation radius is used for this 
dust species instead. 

(3) The disk density distribution, which is either 

(a) described analytically (by the radial density profile and opening angle 
of the disk), or 

(b) provided by the user in form of an upload-file (see Sect. A. 2 for the 
file format). 

(4) The disk mass, i.e. the total mass of all dust grains, 

(5) The dust grain size distribution given by the minimum and maximum 
grain size and an exponent that describes the slope of the distribution, 

(6) The relative abundances of the chemical components. Beside a menu of 
predefined chemical components the possibility to upload further compo- 
nents is provided (see Sect. A. 3 for the file structure). These abundances 
represent either relative mass densities per volume element or relative 
number densities of dust grains. 

(7) The specification of the observable SED to be simulated (wavelength 
distribution). 

B: 

The handing over of the input parameters and file upload are managed by 
PERL and CGI scripts. The input parameters are stored in a file which is 
processed by unit C (see below). 

C: 

The OTC is a radiative transfer tool optimized for optically thin configurations. 
It is written in Fortran 90 and solves following tasks: 

(1) Evaluation of the input data: 

(a) Data type verification for each input value, 

(b) Test of the model integrity, and 

(c) Test of specific upper / lower limits of input parameter values. 

(2) Calculation of the SED: 

(a) Definition of the wavelength range and wavelengths at which the 
dust absorption (and based on this the dust temperature distribu- 
tion) shall be calculated. In the case of a blackbody radiative source, 
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this wavelength range is derived under consideration of the black- 
body temperature. A fixed number of wavelengths (~ 500) is then 
distributed logarithmically equidistantly within this range. In the 
case of uploaded SEDs for the radiative source or the choice of one 
of the SEDs from the internal database (see unit D below), the given 
wavelengths / wavelength range are used. 

(b) Calculation of the dust absorption efficiencies at the wavelengths 
of stellar emission (for dust absorption only) and the user-defined 
observing wavelengths for each grain size (o$) and chemical compo- 
sition 

(c) Calculation of the radial distances r(T, ai,£j) of the dust grains to 
the heating source that correspond to predefined temperatures (see 
Eq. 5; 2.73K < T < sublimation temperature). 

(d) If the density distribution is provided on a grid, the array r(T, a, £) is 
interpolated in order to provide the identical radial spatial resolution 
as the uploaded density distribution. 

(e) If required: Calculation of the scattering efficiency of the dust grains 
at the observing wavelengths. 

(f) Calculation of the relative net contribution of each individual dust 
species outside the corresponding sublimation radius. 

(g) Summation over the net contributions and normalization of the total 
flux based on the total dust mass in the model. 

(3) Creation of all output files: 

(a) for the user: final results (resulting SED, radial temperature distri- 
bution, etc.), 

(b) internally: intermediate results, such as the relative contributions of 
each dust species to the net SED; these files are required in unit F 
(see below). 

D: 

A selection of astrophysically relevant dust data and SEDs of radiative sources 1 
is compiled in an internal database. The optical data of following dust species 
have been included so far: 

(1) Silicates, oxides, and carbon configurations published by Dorschner et 
al. 1995 and Jager et al. 1998 made available at 

http : / /www . astro . uni- j ena . de/Laborat ory/Database/ odata . html 

(Henning et al. 1999), and 

(2) "Astronomical silicate" and graphite published by Weingartner & 
Draine (2001). 



1 So far, only the solar SED is included (from measurements published by Labs & 
Neckel 1968). Further SEDs are planned to be included according to the response 
of the user community. 
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Detailed references are given in unit A. 
E: 

Presentation of all results on a WWW page (see Fig. 1 for an example). 
F: 

The WWW page with the results also allows to modify selected model param- 
eters, such as 

(1) Disk mass, 

(2) Grain size distribution slope, 

(3) Relative abundances of individual chemical components, and 

(4) Distance. 

Based on the intermediate results stored in unit C (step 3), a simple new 
weighting of the relative flux contributions from the individual grains species 
allows a quick calculation of the SED of the modified model. 
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Simulation protocol 



Wavelenqih Imicronl 



| Modify Model | 



Define new Model 



Some remarks on your model 

,., , j 

- Wavelength range of stellar emission (for dust heating): 1.000000E-03 - 10 00.000 micron.. 

Inner; lk the du$i . ■■: r- . : 

- Optical depth in the raidplane(@ shortest wavelength = 1.000000E-03 micron) - 3.745614E-0S 

- Dust reemwfetf^^ii^^-^^fcl-Oiriinosity - to - stellar luminosity ratio: 1 .603834E-04 



Messages providing an insight 
into further characteristics 
of the model and minor 
problems that occured during 
the calculations 



Additional files for download 

1 Spectral Energy Distribution of the Star 

....... . . . 

■J '. . , ■ ■ . - ■ ' ' ■* ' : - ■ 

■: ■ ■ ; . ; . ■ ■ ■■ ■■: ■ ■ 

i 

■ ■■ ■ ■■ ■■ ■ ■ ■■ . ■ '.. ■ : ■ . ■■ .; ■ .. ■ 



Additional results (download] 



Disk Dust Mass 
(T (3.326254E-05 M [M(Earth)] 
Dust Grain Size Distribution 



C (1.000000E-10 M [M(S-lUI)] 



mia grain radius [micron]: 0.1000000 Distribution: n(radius) ~ radiuses, i = 3.500000 

i i i i. , i ,| ... 

Abundances of chemical components frel. mass densities per volume element) 

Silicates 



Possibility to modify 
selected model parameters 



|1 00.0000 



h Astronomical SiHoaM; 



Distance 140.0000 



Questions? Manual 



Main model definition 



Update SEP | 



Fig. 3. DDS: Example output. 



4 Performance of the code 



□ - 

□ 

0.01 O.'O 1.00 10.00 100.00 '000.00 12 3 4 

Grain size [nicron] Grain size distribution 

Fig. 4. Performance of OTC/DDS: Run-time of the code as a function of grain 
size: [A] Single grain sizes; [B] Grain size distributions (1: 0.01-1/xm, 2: 0.1-10/im, 
3: l-lOO/um, 4: 10-1000/um). For this test the dust reemission has been calculated 
for 100 wavelengths linearly distributed in the range 10-100^m. Otherwise, default 
settings have been used, such as 500 wavelengths to simulate the stellar emission. 
The CPU times have been estimated on an Intel(R) XEON(TM) 2.80 GHz processor 
using the Intel Fortran 90 compiler for Linux. 

The runtime of the code is mainly determined by the grain size (distribution) 
and number of wavelengths (a) at which absorption of the stellar light is 
modelled (typically ~ 500) and (b) at which the observable SED has to be 
calculated (defined by the user). The reason for this behaviour is that the 
dust grain parameters (efficiency factors, albedo - see Sect. 2) are calculated 
in each individual simulation in order to provide maximum flexibility in the 
choice of the grain size distribution. For a fixed model setup with a single grain 
size (defined by minimum grain size = maximum grain size) the run-time of 
the code as a function of grain size is shown in Fig. 4 [A]. For selected grain 
size distributions (defined by: minimum grain size < maximum grain size) the 
run-times are shown in Fig. 4[B]. In the latter case, 32 individual grain sizes 
logarithmically equidistantly distributed between the minimum and maximum 
grain size, are considered for each chemical component. 

The OTC/DDS has been tested by comparing its results (radial temperature 
distribution, reemission and scattered light SED) with those obtained with 
the radiative transfer code MC3D (Wolf 2003). The results agree within the 
numerical accuracy, i.e. both codes produce practically the same results. 
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A Input file formats 

The SED of the illuminating and heating source, the dust density distribution, 
and additional dust species can be uploaded in form of ASCII files. These files 
must have tabular structures which are described in Sect. A.1-A.4. The file 
structures are described in an online help webpage, which is connected via links 
to each of the individual upload sections. Since the particular file structure of 
some of the input files may be changed over time, only the links to the help 
pages are given in those cases. 



A.l Stellar SED 

The file structure for stellar SEDs is identical to internally predefined stellar 
SEDs (see there for examples). It is documented at 

http://aida28.mpia-hd.mpg.de/~swolf/dds/dds-manual. html#stellar_sed-upload 
A. 2 Density Distribution 

If no(r, 9, </)) is the arbitrary (optically thin) density distribution, the radial 
density distribution n(r) required for the simulation of the SED can be derived 
as follows (simple averaging): 



The structure of the file for density upload is documented in Tab. A.l. 




(A.l) 



o o 
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# Header with remarks etc. (optional) 




# ... 




Number radial grid points [integer] 


Radial distance r from the source [AU] 


Corresponding relative dust grain num- 


[Boat] 


ber density n(r) [Boat] 



Table A.l 

File structure required for dust density distributions. 



When providing the set of pairs (r, n(r)), the following guidelines should be 
considered: 

(1) The distance between two subsequent radial points should decrease to- 
wards the heating source in order to allow the code (OTC) to resolve the 
increasing radial temperature gradient. 

(2) The step size between two subsequent radial points has to be chosen 
smaller than the typical size of structures in the density distribution 
in a particular distance r from the star in order to prevent skipping 
of local density enhancements, etc. (for the simulation of the spectral 
energy distribution a linear increase/decrease of the density between two 
subsequent radial points is assumed). 

Further Remarks: 

• The quantity n(r) represents the relative number of grains. The absolute 
number of grains is estimated internally based on the total mass of the dust 
configuration, specific dust density and grain size distribution. 

• A radial grid with a logarithmic equidistant distribution of grid points is a 
good choice in the case of " smooth" density distribution decreasing towards 
the radiative source. In the case of a clumpy structure, however, the grid 
has to be adapted to the size of the clumps in the (radial) regions where 
they are present. The same applies if for instance circumstellar disks with 
gaps (due to planet-disk interaction) are considered. 

• If the dust density distribution results from an n-particle simulation, one 
should subdivide the model space in spherical shells centered on the star 
and estimate the mean relative number density of dust grains in each shell 
(total number of grains Ni / volume of the shell Vi). Let n{ri) and n(r wo ) 
be the number densities at the innermost and outermost shell of a model 
subdivided into 100 shells. Let T\ and r wo be the corresponding mean radii 
of these shells. Then, the file prepared for upload of the density distribution 
might have the following structure (see also Fig. A. 2): 
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# Header with remarks etc. (optional) 


102 - number of subsequent lines 


Inner radius 


n(r\) 


T\ 


n(ri) 




n(r 2 ) 






r99 


n(r 9 9) 


noo 


n(r wo ) 


Outer radius 


n(r wo ) 



Table A.2 

Example upload-file for density distributions resulting from n-particle simulations 
(see also Fig. A.2). 

Image not available 

Fig. A.l. Illustration of the subdivision of the model space in the case of den- 
sity distributions resulting from n-particle simulations. The better the radial 
density distribution is sampled - especially in regions with high density gra- 
dients - the higher is the accuracy of the corresponding SED calculated with 
the DDS (see Sect. A.2). Figure available in the complete article - see: 
http: / / aida28.mpia-hd.mpg.de/~swolf/dds / doc / swolf_dds.pdf 
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A.3 Dust Data 



Header with remarks etc. (optional) 




# ... 






Identifier (e.g. chemical composition) [string] 


Specific dust grain density in units of [g/cm 3 ] [integer] 


Sublimation temperature of the dust grains in units of [K] [Boat] 


Number of Wavelengths 


in the file [integer] 




Wavelength [^m] [Boat] 


Re (refractive index) [Boat] 


im(refractive index) [Boat] 



Table A.3 

File structure required for dust data upload-files. 



A dust species of a particular chemical composition is defined by its spe- 
cific (material) density, its sublimation radius, and its wavelength-dependent 
complex refractive index. A large database of laboratory measurements of as- 
trophysically relevant refractive indices is available at 
http : / /www . astro . uni- j ena . de/Laborat ory/Database/ odata . html 
(Henning et al. 1999). 

The file structure for dust data is identical to those accessible through the 
DDS (see there for examples). It is documented in Tab. A.3. 



A. 4 Observed SED 



The DDS allows to upload observed SEDs in order to overlay them to the 
results of the simulation. The file structure for observed SEDs is identical to 
those created by the DDS in order to allow an upload and overlay of simulated 

SEDs as well (see Sect. B.l). It is documented at http://aida28.mpia-hd.mpg.de/~swolf/dds/dds-manua 
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# Header (model description) 


Wavelength [/xm] [Boat] 


Flux [mJy] [Boat] 


0.0 


0.0 



Table B.l 

Structure of the output file with the calculated SED. 



B Output formats of selected files 

B.l Resulting SED 

The file structure for observed SEDs is identical to those created by the DDS 
(in order to allow an upload of simulated SEDs as well). It is documented in 
Tab. B.l. 



B.2 Radial temperature distribution 



# Header (model description) 


Radial distance to the source r [AU] [Boat] 


Temperature [K] [Boat] 



Table B.2 

Structure of the output file with the radial temperature distribution. 



The structure of the file with the radial temperature distribution is given in 
Tab. B.2. For a dust ensemble consisting of grains with different sizes and 
chemical composition, the radial temperature profile is stored for each single 
grain species. Following algorithm to write out the data is implemented: 

for all chemical compositions 
for all grain sizes 

write brief header containing the current 
chemical composition and grain size 
(starting with a "#" sign) 
for r sub < r < outer radius 

write r, T(r) 
end for 
end for 
end for 
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